function plot_convergence_curve()
% GMCM 2025 优化收敛曲线绘制脚本
% 读取收敛数据并生成美观的可视化图表

fprintf('=== GMCM 2025 收敛曲线分析工具 ===\n');

%% 1. 读取收敛数据
% 优先使用详细日志文件，如果不存在则使用简化数据
detailed_data_file = 'result/Problem1/detailed_convergence.csv';
simple_data_file = 'result/Problem1/convergence_data.csv';

if exist(detailed_data_file, 'file')
    fprintf('使用详细收敛日志: %s\n', detailed_data_file);
    data = readtable(detailed_data_file);
    iterations = data.Iteration;
    scores = data.Score;
    phases = data.Phase;
    times = data.Time_Seconds;
    improvement_pcts = data.Improvement_Pct;
    descriptions = data.Description;
    detailed_mode = true;
elseif exist(simple_data_file, 'file')
    fprintf('使用简化收敛数据: %s\n', simple_data_file);
    data = readtable(simple_data_file);
    iterations = data.Iteration;
    scores = data.Score;
    phases = data.Phase;
    times = data.Time_Seconds;
    descriptions = data.Description;
    improvement_pcts = zeros(size(scores)); % 填充零值
    detailed_mode = false;
else
    error('找不到收敛数据文件: %s 或 %s', detailed_data_file, simple_data_file);
end

fprintf('成功读取 %d 个数据点\n', height(data));

%% 2. 数据分析
initial_score = scores(1);
final_score = scores(end);
phase1_score = scores(strcmp(phases, 'Phase1'));
if isempty(phase1_score)
    phase1_score = scores(2); % 如果没有明确的Phase1标记，使用第二个点
end

total_improvement = (initial_score - final_score) / initial_score * 100;
phase1_improvement = (initial_score - phase1_score) / initial_score * 100;
phase2_improvement = (phase1_score - final_score) / phase1_score * 100;

fprintf('\n=== 收敛分析结果 ===\n');
fprintf('初始分数: %s\n', addCommas(initial_score));
fprintf('阶段一结果: %s (改进: %.2f%%)\n', addCommas(phase1_score), phase1_improvement);
fprintf('最终结果: %s (总改进: %.2f%%)\n', addCommas(final_score), total_improvement);
fprintf('阶段二改进: %.2f%%\n', phase2_improvement);

%% 3. 创建收敛曲线图
if detailed_mode
    figure('Position', [100, 100, 1600, 1000]);
    subplot(2,3,[1,2]);
else
    figure('Position', [100, 100, 1200, 800]);
    subplot(2,2,[1,2]);
end

hold on;

% 绘制连接线
plot(iterations, scores, 'b-', 'LineWidth', 2, 'DisplayName', '优化轨迹');

% 不同阶段用不同颜色标记
initial_idx = strcmp(phases, 'Initial');
phase1_idx = strcmp(phases, 'Phase1');
phase2_idx = strcmp(phases, 'Phase2');

scatter(iterations(initial_idx), scores(initial_idx), 150, 'red', 'filled', 'DisplayName', '初始解');
scatter(iterations(phase1_idx), scores(phase1_idx), 120, 'orange', 'filled', 'DisplayName', '阶段一');
scatter(iterations(phase2_idx), scores(phase2_idx), 100, 'blue', 'filled', 'DisplayName', '阶段二');

% 添加关键点标注
text(iterations(1), scores(1), sprintf('  初始解\n  %s', addCommas(scores(1))), ...
     'FontSize', 10, 'Color', 'red', 'FontWeight', 'bold');
text(iterations(end), scores(end), sprintf('  最终结果\n  %s', addCommas(scores(end))), ...
     'FontSize', 10, 'Color', 'blue', 'FontWeight', 'bold');

% 添加阶段分界线
phase1_end = find(strcmp(phases, 'Phase1'), 1, 'last');
if ~isempty(phase1_end)
    xline(iterations(phase1_end), '--', 'Color', [0.7 0.7 0.7], 'LineWidth', 1.5, 'DisplayName', '阶段分界');
end

title(sprintf('GMCM 2025 优化收敛曲线\n总体改进: %.2f%% (%s → %s)', ...
      total_improvement, addCommas(initial_score), addCommas(final_score)), ...
      'FontSize', 14, 'FontWeight', 'bold');
xlabel('迭代轮次 / 阶段');
ylabel('V_{stay} 内存占用');
legend('Location', 'northeast');
grid on;
hold off;

% 设置Y轴格式
ax = gca;
ax.YAxis.Exponent = 0;
ytickformat('%,.0f');

%% 4. 时间-分数关系图
if detailed_mode
    subplot(2,3,3);
else
    subplot(2,2,3);
end
plot(times/60, scores, 'g-o', 'LineWidth', 2, 'MarkerSize', 6);
title('时间 vs 分数');
xlabel('时间 (分钟)');
ylabel('V_{stay} 分数');
grid on;
ytickformat('%,.0f');

%% 5. 改进速率分析
if detailed_mode
    subplot(2,3,4);
else
    subplot(2,2,4);
end

if detailed_mode && any(improvement_pcts > 0)
    % 使用详细日志中的改进率数据
    improvement_idx = improvement_pcts > 0;
    bar(find(improvement_idx), improvement_pcts(improvement_idx), 'FaceColor', [0.3 0.6 0.9]);
    title('实际改进率 (来自详细日志)');
    xlabel('改进点');
    ylabel('改进率 (%)');
else
    % 计算相邻点的改进率
    improvements = zeros(size(scores));
    for i = 2:length(scores)
        improvements(i) = (scores(i-1) - scores(i)) / scores(i-1) * 100;
    end
    bar(2:length(iterations), improvements(2:end), 'FaceColor', [0.3 0.6 0.9]);
    title('每步改进率');
    xlabel('步骤');
    ylabel('改进率 (%)');
end
grid on;

%% 6. 详细模式的额外可视化
if detailed_mode
    % Phase 1 详细分析
    subplot(2,3,5);
    phase1_data = strcmp(phases, 'Phase1');
    if any(phase1_data)
        plot(times(phase1_data)/60, scores(phase1_data), 'ro-', 'LineWidth', 2, 'MarkerSize', 8);
        title('阶段一收敛细节');
        xlabel('时间 (分钟)');
        ylabel('V_{stay} 分数');
        grid on;
        ytickformat('%,.0f');
    end
    
    % Phase 2 详细分析
    subplot(2,3,6);
    phase2_data = strcmp(phases, 'Phase2');
    if any(phase2_data)
        plot(times(phase2_data)/60, scores(phase2_data), 'bo-', 'LineWidth', 2, 'MarkerSize', 6);
        title('阶段二收敛细节');
        xlabel('时间 (分钟)');
        ylabel('V_{stay} 分数');
        grid on;
        ytickformat('%,.0f');
    end
end

%% 7. 保存图表和统计报告
output_dir = 'result/Problem1';
if ~exist(output_dir, 'dir')
    mkdir(output_dir);
end

% 保存图像
saveas(gcf, fullfile(output_dir, 'Conv_Case0_convergence_analysis.png'));
saveas(gcf, fullfile(output_dir, 'Conv_Case0_convergence_analysis.fig'));

% 生成统计报告
report_file = fullfile(output_dir, 'Conv_Case0_convergence_report.txt');
fid = fopen(report_file, 'w');
fprintf(fid, '=== GMCM 2025 优化收敛分析报告 ===\n');
fprintf(fid, '生成时间: %s\n', datestr(now));
if detailed_mode
    fprintf(fid, '数据来源: 详细收敛日志 (%d 个数据点)\n\n', height(data));
else
    fprintf(fid, '数据来源: 简化收敛数据 (%d 个数据点)\n\n', height(data));
end
fprintf(fid, '=== 总体统计 ===\n');
fprintf(fid, '初始分数: %s\n', addCommas(initial_score));
fprintf(fid, '最终分数: %s\n', addCommas(final_score));
fprintf(fid, '总体改进: %.2f%%\n', total_improvement);
fprintf(fid, '总运行时间: %.2f 分钟\n\n', times(end)/60);
fprintf(fid, '=== 分阶段分析 ===\n');
fprintf(fid, '阶段一改进: %.2f%% (%s → %s)\n', phase1_improvement, ...
        addCommas(initial_score), addCommas(phase1_score));
fprintf(fid, '阶段二改进: %.2f%% (%s → %s)\n', phase2_improvement, ...
        addCommas(phase1_score), addCommas(final_score));
fprintf(fid, '\n=== 关键改进点 ===\n');
for i = 1:length(iterations)
    if detailed_mode && exist('improvement_pcts', 'var')
        fprintf(fid, '步骤 %d: %s (改进率: %.4f%%) - %s\n', ...
                iterations(i), addCommas(scores(i)), improvement_pcts(i), descriptions{i});
    else
        fprintf(fid, '步骤 %d: %s (%s)\n', iterations(i), addCommas(scores(i)), descriptions{i});
    end
end

if detailed_mode
    fprintf(fid, '\n=== 详细分析 ===\n');
    phase1_improvements = sum(strcmp(phases, 'Phase1') & improvement_pcts > 0);
    phase2_improvements = sum(strcmp(phases, 'Phase2') & improvement_pcts > 0);
    fprintf(fid, '阶段一改进次数: %d 次\n', phase1_improvements);
    fprintf(fid, '阶段二改进次数: %d 次\n', phase2_improvements);
    if phase1_improvements > 0
        max_phase1_improvement = max(improvement_pcts(strcmp(phases, 'Phase1')));
        fprintf(fid, '阶段一最大单次改进: %.4f%%\n', max_phase1_improvement);
    end
    if phase2_improvements > 0
        max_phase2_improvement = max(improvement_pcts(strcmp(phases, 'Phase2')));
        fprintf(fid, '阶段二最大单次改进: %.4f%%\n', max_phase2_improvement);
    end
end

fclose(fid);

fprintf('\n=== 输出文件 ===\n');
fprintf('收敛曲线图: %s\n', fullfile(output_dir, 'Conv_Case0_convergence_analysis.png'));
fprintf('MATLAB图文件: %s\n', fullfile(output_dir, 'Conv_Case0_convergence_analysis.fig'));
fprintf('详细报告: %s\n', report_file);
fprintf('\n✅ 收敛曲线分析完成！\n');

end

function str = addCommas(num)
% 为数字添加千位分隔符
str = sprintf('%,.0f', num);
end
